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Analysis of Multiple Cracks in an Infinite Functionally Graded Plate 


N. I. Shbeeb, W. K. Binienda and K. L. Kreider 
University of Akron 
Akron, Ohio 

Abstract 

A general methodology was constructed to develop the fundamental solution for a 
crack embedded in an infinite non-homogeneous material in which the shear modulus 
varies exponentially with the y coordinate (thickness). The fundamental solution was 
used to generate a solution to fully interactive multiple crack problems for stress intensity 
factors and strain energy release rates. Parametric studies were conducted tor two crack 
configurations. The model displayed sensitivity to crack distance, relative angular 
orientation, and to the coefficient of nonhomogeneity. 

Introduction 

One of the disadvantages of composites is the mismatch of the thermal expansion 
coefficients between its constituents. This mismatch produces residual stresses, which 
may initiate debonding, delamination, and micro-cracks. For example, application of 
ceramics as a thermal coating for a metal substrate often produce debonding at the 
interface after a small number of thermo-mechanical load cycles. In order to minimize the 
mismatch between the ceramics and metal a new technology was developed. This 
technology allows fully tailored processing of materials and interfacial zones with 
predetermined continuously varying mechanical properties, that are known as 
Functionally Graded Materials (FGM) (see Asish et. al„ 1997 and Holt et. al., 1993). 


NASA/CR — 1 999-208676 


I 



FGM could be described, as two-phase particulate composites wherein the volume 
fraction of its constituents differs continuously in the thickness direction (see Niino and 
Maeda, 1990; Hirano and Yamada, 1988; Hirano et. al., 1988; and Kawasaki and 
Watanabe R., 1990). This implies that the composition profile could be tailored to give 
the appropriate thermomechanical properties. Their physical properties can be determined 
either experimentally or using higher order theory for FGMs developed by Aboudi, 
Pindera and Arnold ( 1997). 

Delale and Erdogan (1983) solved the crack problem for a nonhomogeneous plane. 
The authors considered the plane elasticity problem, in which the material is isotropic, 
has a constant Poisson’s ratio (v), and the Young’s modulus (E) is of an exponential form 
varying in the x-direction, namely, 

E(x) = E {) e p ' (1) 

where (3 is a non-homogeneity constant and E () the Young modulus of the homogeneous 
material. They found that Poisson’s ratio did not have much effect on the resulting stress 
intensity factors. And that the strain-energy release rate at the crack embedded in the 
portion of the medium with higher stiffness is lower than that corresponding to the crack 
tip in the less stiff side of the material. Hence, the crack will grow in the direction of the 
less stiff material. 

Also, Delale and Erdogan (1988) solved the collinear crack problem for two 
dissimilar homogeneous elastic half-planes bonded to a very thin nonhomogeneous layer. 
The elastic properties of the interfacial material varied continuously between those of the 
two semi-infinite planes. The Airy stress function was used in their formulation of the 
solution in which it was assumed that it is composed of two functions, one is associated 
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with an infinite plane containing the crack on the x-axis, while the second is an uncracked 
strip. Their results showed that if the crack location approaches the less stiff material, the 
strain energy release rate increases. 

It can be noticed that multiple oriented crack problems embedded in a non- 
homogeneous infinite plate have not yet been addressed. Thus the scope of this work will 
deal with the general solution to a single and multiple oriented cracks embedded in a 
nonhomogeneous infinite plate. It is assumed that the FGM has a constant Poisson’s ratio 
and the shear modulus is of an exponential form. The solution is valid for both plane 
stress and plane strain. 


General Problem Formulation 

The solution of the mixed boundary value problems for stress intensity factors or 
strain energy release rates at a crack tip is obtained from the perturbation part of the 
problem, see Figure 1 . Before any particular problem is addressed, the general strategy 
of solution is discussed in this section. 

Assume that there are two states of stresses, one is associated with a local coordinate 
system (x r y,) in an infinite plate, while the other is associated with boundaries of a finite 
plate defined in a structural coordinate system (x-y). The crack lies on the X|-axis, which 
is at an angle 0 from the x-axis. In the case of infinite plate problems only the first state 
of stress exists, but for the general problem the total stresses in the local coordinate 
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where. 


>’,) = sin 2 (0)f7 vi (x, y) + cos 2 (0)<X„.(x, y) - 2sin(0)cos(0)r vv (x, >0 
T vrv,^i’ -Vi) = -sin(0)cos(0)o„(x, y) + sin(0)cos(0)cr vv (x, y) 

+ (cos 2 (0) - sin 2 (0))r iv (x, y) 
x = x, cos(0) - y, sin(0) 
y = x, sin(0) + y, cos(0) 


The stress boundary conditions obtained from the perturbation problem are: 

- p,(x,) = lim a T . (jc,, y, ) 

Vj -»0 11 

-P^S* 0= 1 im< t ,U„y,) 

v .-^0 1 1 


(3) 


where pi(X|) and P 2 (x 0 are the normal and shear tractions of the inner crack surfaces. 
Upon substitution of (2) into (3) the boundary conditions becomes: 

" Pt Ui ) = lim (x, , y, ) + lim a< 2 ’ (x, , y, ) (4) 

y. — >0 11 v. -^0 11 


- p 2 (x, ) = lim T v . (x, , y, ) + lim (x, , y, ) 


v, — >0 


(5) 


It is noticed that the principal part will be produced from the first part of (4) and (5). The 
most general form of the stresses are expressed as: 

1 


U| . V, ) = — J , y } ,t)fj(t)dt 


1 h 

<;,Upy,) = — J K'Ss^y^fSW 

““ t/ 

1 h 

a Z .V, ) = — J K™ (x, , y, , t)fj ( t)dt 


(6) 


1 

T \oS x ^yS = —\ K Si( x ] ,y i ,t)f l (t)dt for a < (x,,r) < b 
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where fj are so called auxiliary functions defined as derivatives with respect to X| of 
displacement jumps along the crack. The kernels are expressed as: 


5 ° (*>’>',>') = j X"\a,y,)e* x '-"da 

Kpix^y^t) = J Xf(a, jc, cos 6 - >-, sin o^w'+y^-'da 


(7) 


The expressions of xj- 1 * and xjp depend on the stress and displacement continuity 

of the problem. If xjj^ do not vanish as |a | approaches infinity then an asymptotic 

analysis is done to separate the singular part from the regular. Consequently equation (7) 
can be integrated numerically. 

As |a | approaches infinity equation (7) becomes: 




n cos 6 )+cxi { a j cos 0 - v, sin 0 — / ) 


da 


Substituting equation (8) into (6), the following is obtained: 

— oo a 


o; 1 ; (*,, y,) = — j f,(t)dt 
1 1 ' 1 2/r J ' J 




j sin v, cost) )+i«( V] cos#- 


y '™ 9 - n da]fj{t)dt 


( 8 ) 


(9) 


( 10 ) 
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Further simplifications can be achieved by taking the limit of equation (9) so that the 
first terms in equations (4) and (5) can be determined: 

lim[-i- J X"'(a)e My ' +lalx '- n da\f i (t)dt] (11) 

v,-*o 2n J J 

-oo a 

Assume that X^" have the following asymptotes: 


asy\ = c'lP + id 0 1 a — > 

u V (12) 

asyl = c 0) - id 0 ' «->-<» 

Splitting equation (12) at oc=0 into two parts and making the change of variable for 
the part from to 0 by letting a=-P, adding and subtracting (12) from (11) and taking 
the limit as yj— > 0, the following is obtained: 


^ JT { [X (a) + X 0) (a) - 2c‘" ]cos(a(r - a-, )) + 

[ X (a) - X ‘“(a) + 2 id" ’ ]/ sin(a(r - a. )) } da + lim[— SIP 0 '] 

1 v, — 2 K 


(13) 


where, X^' (a) is the complex conjugate of Xj^ajand the term denoted by SIP l!> is: 


SIP 0 ' = J 2e~ m ' { c ° 1 cos (a(t - a, )) + d 0) sin(a(/ - a, ))]da ( 14) 

0 

The following integral identities can be used to evaluate equation (14) (Abramowitz 
and Stegun, 1964): 


r - a? 

e w ‘ cos (ma)da = — - 

- n~ + m~ 

r ffj 

e a ‘ sin (ma)da = — - 

J n~ + m“ 


(15) 


Hence, it can be shown that (14) becomes. 


c°'v 

L i j -M 


( t - A, )■ + >7 ( t — A, )' + V’i 


(16) 
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where upon taking the limit of (16) and substituting the result into (13) the following 
expression is obtained: 


— T { [*« ‘ + x u' - 2<4”]cos(or(? - x, )) + 
2k Jo 


1 d 0) 

[X"' - X"' + 2/4"]/ sin(a(r - x t ))}da + -f— 

7T t A 


(17) 


In some problems the integrand X'" (a) does not converge rapidly to zero, 
consequently U (where Xj"i s close t0 zero ) is lar § e ’ thus for computational efficiency an 
additional term of the asymptote is taken which is of the order (a 1 ), namely: 


e (U *!" 

<wvl 1 = ^ + / — 
a a 


■a — » +°° 


(18) 


Equation (17) remains the same if the asymptote (18) is subtracted and added at the same 
time. It becomes: 


Jl) 


2 n 


r { [X "’(a) + X-j’ (a) - 2(4" + — )]cos(a(/ - x, )) + 

JO OC 


. 0 ) 


[X"Aa)-X'"(a) + 2/(4 1) +— )]/sin(a(r-x,))}Ja+ 1 


0 ) 


a 

1 7 24" cos (a(r-x,)) l 7 2#;, 1 ’ sin(a(r - x, )) 

+ — — 2 da — da 

97 7 J TV 9 77 J 


n : t — x. 


(19) 


In { a 2 k 0 


a 


To evaluate the last two terms of equation (19) the following identities are used 
(Abramowitz and Stegun, 1964): 


f — sin(a(r - x, ))da = — T 
l a 2 \t 


k (t - x. ) 


1 

[ — cos(a(r - x, ))da = —Ci(U(t - x, )) 

i a 


( 20 ) 




= -(C 0 + log|(/(t - x,)| + J 


COS P - 1 


dP) 
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where, C 0 is the Euler constant. The following expression replaces equation ( 19), 


2 n 


r~ ) e ( " 

f { [*;'» (a) - X (a) + 2/(</<" + ^-)]/sin(a(/ - x, )) }</a 


+ ~~ £ { [X < 1 1 (a) + X ‘J’ (a) - 2c-:; 1 ] cos(a(/ - x, )) }</a 


2;r 


1 r~ 

+ ^rl { [X u ’ <«> ) ■ + * i” («) - 2(4” + — )]cos(a(r - x, )) \da 
2n JU a 


1 d 0) e 0) (t-x) e a) 

+ --2 — 

K t — X, 2 t — X, 7T 


The first part of (4) and (5) can be expressed as: 


+ 


r J 

a 
b 

Ji 


l»m <,(•*.’ >0 = ~ f “ ' J ~ 1 ~ ’~d t + — j *,7 (x, » Ofj ( t)dt 

v i->° n J t - jc, 27T J 


+ 





it - *,) 

.£? 

■ 2 

|f — JC, | 

7T 


/;(/), 

— — dt + 

“j 

r - 

~ 

2tt J 

u 

£ (,) 

*2j 

(f-jr.) 

e (l) 

2 

r - x, 

- 

K 


Ci(U(t - x, ))]/,(/ )dt 


Ci(U(t- X] WjiOdt 


( 21 ) 


( 22 ) 


(23) 


where, 


4"U .0 = f { [x^(a) - X 4(a) + 2i(<’ +^-)]isin(a(r -x, ))}</a 

+ £ {[X4(a)+x4(a)-2c‘ l, ]cos(a(r-x 1 ))}Ja (24) 

r» ef- 1 ’ 

+ JJ{ [X <" (a) + X (a) - 2(c<" + -£-)]cos(a(r - x, )) } da 

A similar procedure is applied procedure to the second terms in equations (4) and 
(5), expressed in the form shown in (10). The terms denoted by Xj y 2, (a)have the 
following asymptotes: 

asy 3 = (a i2) + + c', 2) + id'?' a — > +~ 

asy4 = (-«4 + ibl 2> )a + 4 21 - i'4 2 ’ a -> — 1 
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Simplifications can be made with the use of the integral identities (Abramowitz and 


Stegun, 1964): 


r , «' — m~ 

ac cos (ma)d(X = — 7-7 

J (o‘+or)‘ 

2 nm 


joce sin (ma)da = 


, i 1 . 2 

(O' + 07' ) 


Hence, the second part of equations (4) and (5) becomes: 

1 r (x 2 sin 2 0 - cos 2 0 (t - x,) 2 ) 
=> lim CT ( ^Cj , y j ) — I [ ~ 5 i ■) n 2 x 2 ’ 2 

vi-*o 11 jt J (x.” sin' 9 + cos' 9{t — x, )' )’ 


+ 


+ 


2b\fx t ( t - jc, ) sin 9 cos 9 
(x 2 sin 2 6 + cos 2 6(t — x, ) 2 ) 2 
c: 2 ) jt, sin 6 + 2d! 2) cos9(t — x. ) 


(x . 2 sin 2 9 + cos 2 9(t - x,) 2 ) 


7—7 ]f j(t)clt 


+ ±]kl] ) (x i ,t)f l {t)dt 
in J 


1 r u\~?(xf sin 2 9 — cos 2 9{t — a:,) 2 ) 

=> lim<»(x„y,) = — J[ 

Vi->o 11 n J (x sin 


+ 


+ 


( jc," sin' 6 > + cos' 9{t - x,)')~ 
Ib^x ,(f - x, ) sin 9 cos 9 
(x 2 sin 2 9 + cos 2 9{t - x, ) 2 ) 2 
d, 2 > x, sin 9 + 2 d\ 2 - cos 9(t - x, ) 


(x , 2 sin 2 9 + cos 2 9(t — x,) 2 ) 

1 h 
2n J 




(26) 


(27) 


(28) 


where. 


k? ) {x i s) = £[[X™{a)-Xg\a) + 2i(ab™ +rfJ 2 ) V <K, ‘ in ®]isin(a(/-x I ))}da 


(29) 


+ { [X-j 2) (a) + X} 2> (a) - 2{cta)~ ] + cf )e' aV| sinfl ]cos(a(r - x, )) }da 

It should be noted that (27) and (28) do not contain any singularity and the 
asymptotic expansion is applied only for computational efficiency. 
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Fundamental Solution 


The formulated equations in the previous section will be used to solve the problem of 
radial multiple cracks in an infinite isotropic FGM as depicted in Figure 2. Since dealing 
with an infinite plate only equations (22), (23) and (24) will be used. But before doing so, 
the fundamental solution of a single crack is required. Konda and Erdogan (1994) solved 
the single crack problem using Navier equations. In this work the same problem will be 
solved using Airy stress function and the shear modulus (i will vary exponentially with 
the global y-axis. 

The shear modulus is defined as follows: 


My) — 

M-Wi) = » 0 e p "+ s " 
6 = ycos(0) 
p = ysin(0) 


(30) 


where y, 5, and (3 are real constants and represent the coefficients of 
nonhomogeneity. 

The Airy stress function U(x i ,y i) are defined as, 

<9 2 U 






, , d 2 U 


T,v i ( x p)'i) = - 


c) 2 U 

dx,dy, 


The stresses and strains are related through: 

du, 1 




dx x 8 


e v,v,(^o Ti) = — - 


e „v, (- v i ’ V| ) ~ » ( 


d}\ 8/i(x,, y, ) 

1 , d w, d v. 1 


[(* + iV v,,, + fa - 3>? Vt V| ] 
i( K - 3K,v, + (*■ + »Kv, 1 


) = 


2 d y, dx x 2pi{x ] ,y ] ) 


(31) 


(32) 
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where K is defined as, 


k = 3 - 4v for plane strain 

K = — — — for plane stress 

1 + v 


(33) 


and the compatibility equation is defined as: 

d 2 e xx d 2 £ yx d 2 e y 

A i-M > l .M 2 - 


d y; dx s d y, 


= 0 


(34) 


From the compatibility equation the following fourth order governing equation for 
U(X|,yi) is obtained, 

V 4 t/ (x, , y, ) - 2d-^—(V 2 U (jt, , _y, )) - 2(5-^— (V 2 U (x, , y, )) + 


J 1 


8 (56 d^Ujx^y^) | p 2 ^ K-3d 2 U(x l ,y } ) ^ d 2 f/(X|,y,) ^ | 

1 + k dx t d y\ 1 + k dy 2 dx 2 


s -, f d 2 U(x x ,y^) K-3d 2 U(x ] ,y i ) 

0( r — ; + - T-r; ) = U 


dy 


l + K d j 


(35) 


where 


v=( & + $ ) 


(36) 


Defining the Fourier transform as follows. 


V(a,y,)= ^U(x x ,y\)e ,w 'dx x 


(37) 


and applying equation (37) to (35), the following characteristic equation is obtained, 


n 4 — 26ir + 


(5 2 Kj — 3 + ^2 _ 2a(i(5 + a) 


IC n + 1 


W + ad 


8() 


v K 2 + 1 


■ / + 2 a ri 


+ a' 


3 - K , 


a 2 + <5 2 1 1 + (5(2 ia - (5) 

K-, + 1 


(38) 


= 0 
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the roots of which are: 


fi\ = — 
2 


Wo = — 

o 


= — 


n, = - 


5 + a 


8-p, 


8 + p, 


8-P 


f 3 — KT 

fC, + 1 

[3 - K 

K + 1 


3 - k 


K + 


' 2)1 

H 


s + p, 


8-P 


3- K 

fC + 1 

3 - K 

K + 1 


V 


+ 4 (a 2 + ia 


P-8, 


+ 4 (or + ia 


p + S 


3 - K 

K + 1 

3 - K ' 

JC + 1 


+ - 
0 


a + 0. 


j 3 - K 

K+ 1 


V 


+ 4 (a 2 + /a 




3 — K 
K+l 


3- v 
*: + 1 


1 

+ - 
2 


8-P, 


3 - K 

K + 1 


V 


+ 4(a 2 + ia 


0 + ^ 


3 — K 

K + 1 


( 39 ) 


The solution to the ODE becomes: 


V(a, v, ) = B, (a)e'" y ' + B 2 (a)e" 2 '' + Bpa)e""' + B 4 {a)e" 4> ' 


(40) 


so that, 


| 00 

U( a-, , .v, ) = — J [B, (o)e ¥l + B 2 (a)c""’ + B,(a)e* >y ' + BJa)e" 4 " ]e'" a da (4 , ) 


Bounded form of equation (4 1 ) can be obtained upon examination of the roots of the 
characteristic equations. The real part of ni and n 2 are negative while that of n 2 and n 4 are 
positive as a approaches infinity. Hence U 2 is defined for positive and negative y, as, 


U 2 (x t , 
U 2 (x ,, 



+ B 2 e" 2 " )e‘" a da ; 
+ B 4 <?" jV| )e ix ' a da\ 


> 0 


y, < 0 


(42) 


Normal and shear stresses must be continuos at yi=0. The continuity conditions can 
be represented by equation (3 1 ) as: 

dU( jc, ,0+) _ dU(x lt 0) 

d >\ d >\ (43) 

U ( a'| ,0 + ) = U{x v 0‘) 
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where, 0 + is for y |>0 and 0 is for yi<0. Using conditions (43) we can eliminate B 4 


and B 4 : 


B ^n±zJtL B]+ n±z^ B ^ 

n 4 — n } n 4 ~>h 
n 4 - n 3 n 4 - n , 


(44) 


The remaining two unknowns can be expressed in terms of the auxiliary functions: 


f t (x ,) = — — [«,(x,,0 + ) - w, (*,,()“)] 

° x \ 

~\ 

/-> (jf, ) = - — [v, U, ,0 + ) - v, (jc, ,0” )] 
d x. 


(45) 


The final expressions of the stresses are obtained using (31) and Hooke’s law, 


<T v,, l Up.v, + ) = T-j[« 


2k 

ll ^ 


2 F { {a)h v - F 2 (a)h ]2 

h \ | /? 2 ■> ^ 12^21 


F\ (a)h 2] + F 2 (a)h u ^„ 2V| 


K ] /?22 2 ^2 1 


(46) 


^ 1| v 1 (^<>', + ) = -^ j« 2 c 


1 7 2r F } (a)h 22 - F 2 (a)h n 


, «i v, 


//,,/»„ - /i n /i. 


+ 


F t (oc)h 2 \ + F 2 (a)h ] [ ^n,\| 


h u h 22 - h ]2 h 2 , 


(47) 


T,,(^y;> = -TT J«[» 


F t (a)h 22 - F 2 {a)h n 


,»1>| 


h u h 22 - h n h 2l 


+ n-, 


2k 

- FAa)h 2] + F 2 (a)h 


h) \ F 22 h\2^2\ 


e n '']e‘ x ' a da 


(48) 
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where. 


, K+ 1 

«n =— — (». -«0(« i -« 4 ) 

»/i 0 

K* + 1 

/i p =— — (n, -«,)(«, — n 4 ) 
8/i 0 " 


/?,, = 


icx-B , («. -«,)(«. — «, ) 

— — * -(g-(l + y)-g-(ir 2 -3))[ — ' —r — ' — ] 


8/4 


= ^(cr( I + v) - SHk, - 

8/i ( , (<S -« 2 )(5 -/i,)(<5 -n 4 ) 

b 

Fj(a) = j fj(t)e iP ~ ia) 'dt j = 1,2. 


( 49 ) 


The singular integral equations can be solved for the auxiliary functions using the 
boundary conditions: 


- />,(*,) = lim tr (X| , y, ) a < x, < b 


V| — >0 

- p 2 ( jc, ) = lim T (jfp y, ) a < x, < b 

v , — >0 11 


(50) 


Equation (47) and (48) are rearranged as follows. 


l r -.fi e" |V| — h e"''' 

^„„(-v-.vr) = -— J«-[ • - - ■ f, ( «) 


\\ ll 22 2^21 


+ 


Ke™' -h„e 


"|.V| 


h h — h h 
rt \ \ ri 22 ri \ 2 ri 2 \ 


F 2 (a)]e ix ' a da 


(51) 




n \h 22 e n °' — n 2 h 2] e" 
h\ | h 22 /i | ? /i 2 ! 


F,(a) 


+ 


n 2 h u e" 2 '' - nji n e"' y ' 
b\\b 2 2 b\ 2 h 2 \ 


F 2 {a)]e"' ,x da 


(52) 


As a goes to infinity ni=n 2 so Xjj can be expressed as defined in equation (9): 


*„ = 


* 2 , = 


-cr[ 

-ia[ 


h-i-t ft-. 


/i| \h 22 /i] t ^2 1 

/i|/i 22 - n 2 h 2l 

fl\ j /l , 2 ^12^21 


] = -« 2 [ 


^11 ( ? 12 


/i| 1^22 


/i|2^2l 


] = -/«[ 


, ~ Mi/ip 

/i| 1^22 /i 1 2 /i 2 1 


(53) 
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the asymptotes of Xjj are found as: 


5 S 

X, ,(+-) = X,,(-oo) = 

a(l + k) - a(l + k) 

2 i B .. . - 2 / 


X,,(+°°) - 


X ->,(+<») = 


+ 


(1 + k) a{\ + k) 


2 1 


+ 


P 


Xp(-oo) = 


X,.(-oo) = 


p 


(1 + k) - a( 1 + k) 


- 2 / 


P 


( 54 ) 


(1 + jc) a(l + v) 

-S -8 

X, 1 (+°°) - — 7 X 22 (— °°) — 


(1 + k) - a( 1 + k) 


a(\ + K) -a(l + K-) 

that from (54) it is concluded that d ij’ e ij ,C ij and g ij described earlier in equations 
(12) and (18) become, 


so 


= 4 = 

2 

1 + K 

1 

6 

— ^ 22 " 

1 + K 

- J - 

P 

— ^21 “ 

1 + K 

1 

~ C \2 ~ ' 

l 

c 21 = C 


(55) 


Now substituting (53) and (55) into (22) and (23) the final singular integral 
equations are reached, 

K + 1 , . If / 2 (0 


— p,u,) = - f + - f k r ( Xl ,t)MOdt 

2juUpO) * J f - n{~ 

- - - x x ))]f 2 (t)dt 

^ a 

h 


(56) 


i b <> 

+ - f ft,, ( jc, , t) {t)dt - - f ^Ci(U(t - jc, ))]/, {t)dt 

K : 2 

p , (jr ) = 1 \MLdt + - f ft„(jc„ t)Mt)dt 

2 jU (.v„0) / - 7T J / - x, jrJ ' 

- - l^CKUit - x f ))]Mt)dt 


n J 2 


( 57 ) 


i b o 

+ — [ k„(x t ,t)fi(t)dt + — \ —Ci(U(t - x x ))]f,(t)dt 
;r J K J 2 
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where. 


*ii(*i.O = “ X,,]/sin(a(r - x,))}</a 

+ J 0 {[X M + X ll( ]cos(a(r - x^)}da (58) 

+ ru^n + x tu -2— -]cos(a(r - x s )))da) 

Jv «(1 + K ) 

A- |7 (.v,,/) = - + - } [” {[X P( . -X p + — ]/sin(a(/ - x,))jda 

4 jo - 1 + k 

+ j[ {[X 12 + X l2l ]cos(a(f - x t ))}da (59) 

+ Jj[ x n + X Ut -2 — — ^ — -]cos(a(r - x t ))}da] 

Jo a(l + k ) 

A : ,(a-|,/) = - j 1 {[ {[X 21( . - X,, + 4< ]/sin(a(t - x ] ))}da 

4 jo 1 + k 

+ | ( {[X 2I + X 2I< ]cos(a(r - x^da (60) 

+ f { [X,, + X,, - 2 ] cos(a(; - x. )) )da) 

Ju Cf(l + k) 

k 22 (x ] ,t) = — + 1 ( f {[X 22( . - X 22 ]/ sin(a(r - x A )))da 

^ J() 

+ f {[X 22 + X 22 Jcos(a(r -*,))}</« (61) 

JO 

OO 5J 

+ [ { [X„ + X, 7 + 2 ] cos(a(r - jc. )) } da] 

J,! " " a{ 1 + k) 

The definitions for the stress intensity factors (SIF) and the strain energy release rate 
(SERR) can be found in Konda and Erdogan (1994). Applying the Lobatto-Chebyshev 
collocation integration technique, as in Binienda and Arnold (1995), to the system of 
singular integral equations (56) and (57), the normalized mode I SIF were produced and 
compared to that of Konda and Erdogan (1994) as summarized in Table I. 
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Multiple Crack Formulation 

To formulate the multiple crack problem the total stresses of the system needs to be 
determined. The cracks are located along their local Xi axes, which are related by the 


following coordinate transformation: 


A', = x M cos(0, +l - 0,) - y, +1 sin(0, +1 - 0,) 
y, = X i+l sin(0, +1 - 0,) + V, +I cos(0, +1 - 0, ) 

x, +l = x, cos (0, +l - 0,) + y, sin(0, +l - 6 i ) 

y, +1 = -a-, sin(0, +l - 9, ) + y, cos (0, +l - 9,) 

0 M > o, 


(62) 


and the stresses are related through the Cauchy stress transformation tensor: 



o, t , ' 


m 2 
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n ~ 

- 2 mn 
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2 mn 
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0, . 


mn 
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-> o 

m~ - n 
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k..J 
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m " 
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<7,, 
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I) , | V )t ] 

> = 
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n " 
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- 2 mn 

< 

ff.w, 

► 

T 

• v i+i v (>i 


- mn 

mn 

7 7 

wT — ti 


0-v,v, , 



) L -J k J 

m = cos(0, +l - 0,);w = sin(0 /+1 - 0,) 


The material constants are: 


A = 7 sin(0, ) 
<5, = y cos(0 ; ) 


The stresses for each coordinate system are expressed as: 


<7 ,,,u,,.vr)=^- j[ 
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(64) 


(65) 
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( 66 ) 
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where. 


( ' ) is- i j (/) 0) O) (/) 

N u =— (n,-n,)(n 3 -n 2 ) 

(/) K* + 1 (/) <') 

/V |2 = («4 _ ”2K ,J l _n 4) 


0) 

*2.= 


</> 

(V-n = 


8^o 

ia -P,„2 


8jU 


(cr(l + fC 2 )-<S,-(?e-3))[ 


(/) (0 (0 <n 

(n,-n 3 )(n 3 -tt 2 ) 

(I) 0) O') 

(5, -«,)( 5, -«,)(<5, -n 3 ) 

0) O') 0) </) 


] 


8^o 


, o,, , ov , r (/J 4 -nO(Ar,-« 4 ) n 

(a-(l + fe 2 )-5,-(K--3))[ — j77 rrr* ] 


0 > 


(5, -n 2 )(d, -/!,)(< 5, -« 4 ) 


( 71 ) 


and hjj are as defined in (49) for the i ,h crack. Assume that there are n cracks present, then 
the stresses for the i th crack could be expressed as: 


n 

= cr; v (x,,.y,) + ^(T^Ix/a-,,^), v^x,, v,)J 




,<v, (A > >'. ) = T >,v,Uo >/) + S T V,V, [*;(*,’ >0’ V/A’ >’/)] 


(72) 


>=<+i 


u> (./) 

where v, » T >,>•, are found using (63) and they are evaluated as in (9) and (10) 
respectively. One must be carefull when chossing the stress components for the positive 
yi (yi + ) and for the negative y t (yf). Thus the final singular integral equation could be 
expressesd for the i lh crack as: 
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«-+ 1 (i) i > f it \ if *<■> <'> 

, (-TT^ pM)) = L 7^', + - J*2, u, 
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(74) 


<;> (/) u) _ _ << 

where the constants a •••</ and & are defined in the Appendix, the kernels 


are defined as in (58) through (64) and 0=0j-0j. 
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The above solution is reduced to the case of two collinear cracks embedded in the 


isotropic plate to demonstrate high accuracy of the results as shown in Table II. 

Parametric Studies 

In the following parametric studies the length of all cracks is chosen to be 2c = 2. 
The infinite plate is subjected to normal stress along y global direction, a y y =l psi. Cracks 
are located along their local Xj axes, which can be inclined with respect to the global x 
axis. All the geometrical dimensions are normalized with respect to c. The parametric 
studies are presented for the normalized mode-I and mode-II SIF, i.e., k|/k ( , and k 2 /ko, and 

8 k ~ 

normalized SERR, i.e., G|/Go and G-./G () , where k„ = c vv Vc and G„ = — 0 . 

y> 7t( K + 1 ) 

First study takes into consideration the problem of two collinear horizontal cracks. 
The same crack configuration was used to produce the results in Table II for 
homogeneous materials, see insert in Figure 3 or 4. Here we extended the material 
properties to the show transition from homogeneous material to FGM. The distance 
between the inner tips is denoted by r. Figures 3 and 4 show mode I normalized SIF and 
mode-I normalized SERR versus a normalized crack tip distance r/c. The curves are 
shown using the logarithmic scale for the crack distance variable. 

It can be noticed that as the distance between the cracks becomes smaller the SIF 
and SERR become larger for every power of the exponential variation coefficient y. 
Note, that the homogeneous case is represented by making the coefficient y = 0. Both 
driving forces increase as y increases for materials becoming more nonhomogeneous. The 
increase is especially significant for the crack tip distance less than 0.01. 
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The case of collinear inclined cracks at 30 degrees from the horizontal axis is 
shown in Figures 5 to 8. Mode-I SIF shown in Figure 5 has larger magnitudes for both 
the homogeneous and FGM materials than the corresponding Mode-II SIF shown in 
Figure 6. Both modes show increase of SIF for decreasing of the crack tip distance and 
for increasing of the coefficient y. Figures 7 and 8 display mode-I and II SERR. 
Magnitudes of Mode-I SERR are almost three times larger than the corresponding mode- 
II SERR. Both SERR modes increase for decreasing of the crack tip distance r (same as 
SIF) and for decreasing of y (opposite to SIF). It should be pointed out that the material 
stiffness at each crack tip dominates the results for SERR to the point of reversing the 
trend in comparison with SIF with respect to y. 

The cases when two cracks are located along two different local radial axes 
distance d = 1 from the origin is shown in the remaining figures. The location of the first 
crack is kept constant at 30 degrees while orientation of the second crack is changed from 
45 to 90 degrees. Both SIF and SERR are shown for each crack tip versus the orientation 
angle of the second crack. 

Figures 9 and 10 display mode-I and II SIF while Figures 11 and 12 represent 
both modes of SERR at the left crack tip of the stationary crack 1. It can be noticed that 
when crack 2 comes closer to crack 1, the tip ai is shielded and all driving force 
components are significantly reduced. Mode-II SIF has its maximum for the orientation 
angle of the second crack close to 70 degrees. By increasing y higher magnitudes for 
mode-I SIF and lower magnitudes for mode-II SIF are produced. 

Both modes of the SERR depend not only on square of the SIF but also on the 
material stiffness at the crack tip. This influence is especially shown in Figure 1 1 where 
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the homogeneous case produces the smallest SERR for angle 45 degrees similar to the 
corresponding SIF. However, homogeneous SERR curve becomes largest for higher 
angles, which is opposite to the corresponding SIF. 

The opposite character of SIF and SERR is even better shown in Figures 13 to 16. 
Here the shielding effect does not exhibit itself. Both modes of SIF for the homogeneous 
material are smallest (see Figures 13 and 14 for mode-I and II SIF), while both 
components of SERR are largest (see Figures 15 and 16 for mode-I and II SERR), 
because of the crack tip material stiffness influence. 

The results for the left crack tip of the second crack are shown in Figures 1 7 to 20. 
Mode-I SIF depends on the crack orientation and for homogeneous case is smallest. The 
negative values of the mode-I SIF should be interpreted as the crack closure and SERR 
for such case is zero. 

Mode-II SIF at the tip of the second crack versus its orientation is shown in 
Figure 18. For the homogeneous case the maximum k 2 is at 45 degrees. In the cases of 
higher coefficient y the maximum k 2 is shifted towards 60 degrees crack orientation. The 
magnitudes of SIF are larger for nonhomogeneous cases than for the homogeneous case 
when the orientation angles of the second crack are higher than 60 degree. 

Mode-I of SERR monotonically decreases to almost zero at 75 degrees and at 
about 50 degrees does not depend on y, see Figure 19. For 45-50 degrees the 
homogeneous case is the highest and for 50-90 is the smallest. Mode-II SERR is the 
highest for 45-90 degrees angle orientation and has its maximum at 45 degrees, see 
Figure 20. The shift of the maximum to 55 degrees can only be observed for y = 1.0. 
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Mode-I SIF for the right crack tip bi is very small for all the angles examined or it 
is negative (crack closure) for the angle more than 65 degrees for the homogeneous case 
and down to 50 degrees for FGM with y=l. Mode-II SIF is shown in Figure 22. The 
maximum is moved to 60 degrees because of the influence of the crack below. The 
homogeneous material case is the smallest for all the angles of the crack orientation. 

Mode-I SERR at the same crack tip is shown in Figure 23. Homogeneous case 
starts to be the highest at the orientation of 45 degrees and quickly goes to zero at 62 
degrees. The FGM with the highest y also goes to zero but at the smaller angle because 
the crack tip remains closed. Mode-II SERR is shown in Figure 24. All the curves have 
their maximum at 55 degrees. The highest SERR is for y= 1.0 and the smallest for y = 
0.25 at the angle of 45 degrees. The homogeneous material produces the smallest SERR 
for the orientation close to 75 degrees. 

All of the above parametric studies demonstrate that the effect of the material 
properties, crack orientation, location of the additional crack are interdependent and 
consequently produce behavior different than isotropic homogeneous materials. The 
model developed in this work can be used to study fracture problems in FGM and can be 
used to tailor the properties in order to reduce driving force components and effectively 
increase live of these materials. 

Conclusions 

Application of the general solution to the mixed boundary value problem was 
demonstrated to provide an elegant way of obtaining the fundamental solution for a crack 
embedded in an infinite nonhomogeneous plate. The fundamental solution was used to 
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address multiple crack problem. Parametric studies for the multiple cracks revealed, that 
both SEF and SERR highly depend on the crack geometrical parameters such as crack 
orientation, location, relative distance, etc., but they also depend on the power of the 
exponent describing the rate of change of the material elastic parameters, y, and local 
stiffness of the material at each crack tip. 

The results demonstrated that the driving forces can be amplified by the collinear 
crack orientation or they can be reduced by the shielding effect between cracks above or 
below. The character of the amplification or shielding remains similar for 
nonhomogeneous materials but in most cases higher than zero coefficient y increases SIF 
and reduces SERR. 

The well known one-to-one relation between SIF and SERR curves is not always 
valid for FGM, because SERR also depends on the material elastic constants. Hence, the 
SIF curves may have different character than SERR curves. The application of the 
driving forces to crack propagation criterion need to be further studied to determine 
which driving force (SIF or SERR) best correlates with appropriate experimental results. 
However, since SERR includes the influence of SIF and material stiffness at the tip it is 
recommended that total SERR is used for the driving force parameter for FGM. 
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Table I Verification of the Solution. 


cy 

Konda and 
Erdogan (1994) 
k|(a)/Vc 

Present Study 
k|(a)/ Vc 

Konda and 
Erdogan ( 1994) 
k 2 (a)/ Vc 

Present Study 
k 2 (a)/ Vc 

0.25 

1.036 

1.036 

0.065 

0.062 

0.50 

1.101 

1.101 

0.129 

0.122 

1.0 

1.258 


0.263 



c=(b-a)/2 


Table II Two Collinear Cracks In Isotropic Plate 



From literature 



r j_ 

Horri and Nemat-Nasser 

Erdogan 

Present Method 

c 

(1985) 

(1962) 




Inner 

Outer 

Inner 

Outer 

Inner 

Outer 

0.22 

— 

— 

1.45387 

1.11741 

1.45736 

1.11786 

0.50 

1.2289 

1.0811 

1.22894 

1.08107 

1.22894 

1.08107 

0.857 

1.1333 

1.0579 

1.13326 

1.05786 

1.13329 

1.05787 
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Figure 1 . Methodology of Solution for the Fundamental Problem, 
(a) Mixed Boundary Value Problem for the FGM. 

(b) Infinite FGM Plate without Crack. 

(c) Perturbation Problem of a Crack loaded by Surface Tractions. 
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Figure 2. Multiple Cracks Embedded in the Infinite FGM Plate. 



Figure 3. Mode I normalized SIF versus normalized inner crack tips distance 
for two collinear horizontal cracks. (G yy = 1 .0, G xx = o xy =0.0). 
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Figure 4. Mode I normalized SERR versus normalized inner crack tips distance 
for two collinear horizontal cracks. (a yy = 1.0, a xx = a xy =0.0). 
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Figure 5. Mode I normalized SIF versus normalized inner crack tips distance 
for two collinear cracks along the 0=30 deg. (a yy = 1.0, g xx = o xy =0.0). 
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Figure 6. Mode II normalized SIF versus normalized inner crack tips distance 
for two collinear cracks along the 0=30 deg. (o yy = 1.0, a xx = O xy =0.0). 
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Figure 7. Mode I normalized SERR versus normalized inner crack tips distance 
for two collinear cracks along the 0=30 deg. (G yy = 1 .0, a xx = a xy =0.0). 
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Figure 10. Mode II normalized SIF at the tip a, versus crack (2) orientation angle 0 2 for 
two inclined cracks. (o yy = 1.0, o xx = a xy =0.0, d=1.0, 0i=3O deg.). 



Figure 1 1. Mode I normalized SERR at the tip a, versus crack (2) orientation angle 0 2 
for two inclined cracks. (a yy = 1.0, o xx = a xy =0.0, d=1.0, 0i=3O deg.). 
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Figure 12. Mode II normalized SERR at the tip versus crack (2) orientation angle 0 2 
for two inclined cracks. (a yy = 1.0, a xx = a xy =0.0, d=1.0, 0i=3O deg.). 



Figure 13. Mode I normalized SIF at the tip b| versus crack (2) orientation angle 0 2 for 
two inclined cracks. (a yy = 1.0, a xx = a xy =0.0, d=1.0, 0|=3O deg.). 
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Figure 14. Mode II normalized SIF at the tip b| versus crack (2) orientation angle 0 
for two inclined cracks. (o yy = 1.0, o X x = cj xy =0.0, d=1.0, 0|=3O deg.). 
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Figure 15. Mode I normalized SERR at the tip b] versus crack (2) orientation angle 0 
for two inclined cracks. (o yy = 1 .0, G xx = o xy =0.0, d=1.0, 0|=3O deg.). 
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Figure 16. Mode II normalized SERR at the tip b| versus crack (2) orientation angle 0 2 


for two inclined cracks. (c y y = 1.0, a xx = a xy =0.0, d=1.0, 9|=30 deg.). 
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Figure 17. Mode I normalized SIF at the tip a 2 versus crack (2) orientation angle 0 2 
for two inclined cracks. (c yy = 1.0, a xx = a xy =0.0, d=1.0, 0i=3O deg.). 
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Figure 22. Mode II normalized SIF at the tip b 2 versus crack (2) orientation angle 0 2 
for two inclined cracks. (a yy = 1.0, a xx = <5 xy =0.0, d=1.0, 0i=3O deg.). 



Figure 23. Mode I normalized SERR at the tip b 2 versus crack (2) orientation angle 0 2 
for two inclined cracks. (o yy = 1.0, a xx = a xy =0.0, d=1.0, 0|=3O deg.). 


N AS A/CR— 1999-208676 


47 



Figure 24. Mode II normalized SERR at the tip b 2 versus crack (2) orientation angle 0 
for two inclined cracks. (o yy = 1.0, g xx = a xy =0.0, d=1.0, 0i=3O deg.). 
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